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Figure  1.  Single  Doppler  velocity  pattern  of  a  distant  and  stationary 

cyclonic  vortex  as  seen  by  a  Doppler  radar  scanning  with  Infinite 
resolution  at  low  elevation  angle.  2 

Figure  2.  Block  diagram  of  pattern  recognition  algorithm.  5 

Figure  3.  Decision  boundary  in  the  shear-momentum  space.  The  thresholds 

are  indicated.  6 

Figure  4.  Decision  boundary  in  the  velocity  difference  -  azimuthal  distance 

space  with  actual  thresholds  used  by  the  algorithm.  7 

Figure  5.  Solidly  rotating  mesocyclonic  core.  Quantities  used  in  the  text 

to  define  various  calculations  of  diameters  are  indicated.  9 

Figure  6.  Normalized  Doppler  velocities  of  a  combined  Rankine  model  vortex 
centered  at  a  radar  range  r.  Vortex  radius  Is  r+  and  the  ranqe 
resolution  is  much  finer  than  the  size  of  the  mesocyclone. 


11 

Figure  7.  Azimuth  versus  range  display  of  mesocyclone  position  and  diameter. 

Bars  with  circles  depict  the  diameter  D-|  and  position  produced  by 
the  algorithm,  whereas  bars  with  squares  are  values  subjectively 
calculated  from  Doppler  velocity  fields.  Maxi-tornado  damage  path  is 
stippled  area.  13 

Figure  8.  Plot  of  objectively  determined  average  shear  and  subjectively  deter¬ 
mined  maximum  shear.  13 

Figure  9.  Mesocyclone  feature  detected  by  the  algorithm.  Velocities  are 

* 

printed  in  an  azimuth  (on  abscissa)  versus  range  (ordinate)  format. 
Elevation  angle  is  1.2°.  Solid  lines  bracket  all  detected  vectors 
of  a  feature.  Dashed  lines  correspond  to  two  nondetected  vectors.  14 
Figure  10.  Same  as  Figure  9  but  at  an  elevation  of  3.3*.  Two  separate  features 
shown  here  were  detected  in  this  scan.  15 

Figure  11a  The  position  and  sizes  of  features  detected  at  two  consecutive 
scans  and  same  elevation  3.3°.  (1  and  2)  separated  by  6  minutes. 

Azimuthal  diameters  are  indicated  with  vertical  lines  and  range 
diameters  with  horizontal  lines.  Two  features  A  and  o  were 
detected  in  both  scans.  The  two  features  in  the  first  scan  (1) 
correspond  to  the  ones  on  Figure  9.  17 

iv 


Figure  71b 
Figure  77c 
Figure  A.l 
Figure  A. 2 
Figure  A. 3 
Figure  A. 4 
Figure  A. 5 
Figure  A. 6 
Figure  A. 7 


Same  as  77a  but  at  the  next  elevation  angle  (4.5°). 
Same  as  17a  but  at  5.4°  in  elevation. 

Flow  chart  of  the  Main  Program. 

Flow  chart  of  the  Subroutine  SELECT. 

Flow  chart  of  the  Subroutine  FEAT. 

Flow  chart  of  the  Subroutine  CALC. 

Flow  chart  of  the  Subroutine  RDREC. 

Flow  chart  of  the  Subroutine  PRELIM. 

Flow  chart  of  the  Subroutine  UNFOLD. 


1.  INTRODUCTION 


Doppler  radar  is  capable  of  making  accurate  estimates  of  the  radial  compon¬ 
ent  of  velocities  within  severe  thunderstorms.  A  trained  observer  can  examine  a 
field  of  radial  velocities  and  recognize  the  characteristic  signature  of  various 
flow  patterns.  Donaldson  (1970)  stipulated  criteria  whereby  a  vortex  can  be 
identified  from  single  Doppler-radar  observations.  Briefly,  there  must  be  a 
localized  region  of  persistently  high  (>5x10  s  )  azimuthal  shear  (i.e.,  the 
velocity  gradient  along  an  arc  at  constant  range)  which  has  a  vertical  extent 
equal  or  longer  than  its  diameter.  These  criteria  were  further  refined  by 
Burgess  et  al .  (1979)  who  used  them  successfully  to  detect  mesocyclones  and 
inform  the  National  Weather  Service  arid  the  USAF  Air  Weather  Service  of  severe 
thunderstorm  potential.  The  manual  process  consists  of  recognizing  circulation 
patterns,  confirming  height  and  time  continuity  of  the  patterns,  and  verifying 
that  the  radial  velocity  shear  across  the  circulation  is  larger  than  a  range 
dependent  threshold  value.  This  report  concentrates  on  the  automatic  recognition 
of  circulation  patterns  from  Doppler  velocity  fields.  Confirmation  of  height 
and  time  continuity  is  not  considered;  however,  identification  of  large  shear 
is  an  inherent  part  of  the  technique  and  is  discussed  in  detail. 

2.  THE  MES0CYCL0NE 

It  has  been  known  for  some  time  that  a  cyclonic  circulation  larger  than  the 
tornado  may  exist  in  thunderstorms  (Brooks,  1949).  The  larger  rotation  precedes 
the  tornado  and  acts  as  a  vorticity  producing  source.  Fujita  (1963)  labeled 
this  structure  the  mesocyclone,  noting  that  its  diameter  is  too  small  to  be 
observed  with  conventional  surface  weather  observations .  Previous  Doppler  studies 
have  shown  that  single  Doppler  radar  can  easily  detect  and  track  mesocyclones 
(Donaldson,  1970;  and  Burgess,  1976).  From  Doppler  data,  the  average  mesocyclone 
is  known  to  last  for  over  an  hour  with  a  rotational  velocity  of  23  m-s  a 
core  diameter  of  5  km;  and  a  vertical  extent  of  8  km  [further  characteristics 
are  contained  in  Burgess  et  al .  (1982)]. 

The  flow  in  a  mesocyclone  theoretically  can  be  modeled  as  a  Rankine  combined 
vortex,  characterized  by  two  flow  regimes.  The  Doppler  velocity  field  mapped  by 
an  idealized  Doppler  radar  (perfect  resolution)  is  depicted  in  Figure  1.  The 
inner  region  is  the  vortex  core  where  velocity  is  proportional  to  radius  at  all 
radii  out  to  the  radius  of  maximum  tangential  velocity  (solid  body  rotation). 

The  outer  regime,  outside  of  the  core,  features  velocities  which  decrease  with 
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" - °'4  3.  PATTERN  RECOGNITION  TECHNIQUES 
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v  Many  schemes  for  recognizing  two 

dimensional  patterns  have  been  utilized 
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in  various  fields  of  science.  These 
,  ..  ,  can  be  grouped  into  the  following  three 


T?  * _ _  7  r*  •  7  rv  1  -7  ■  ,  ,  wuri  yi  uupuvi  1  1 1  uv  ui  ic  lull  \Jrt  mu  uii  cc 

F-igwe  1.  Single  Doppler*  veloctty  pat -  r  3 

tern  of  a  distant  and  stationary  cyclonic  distinct  categories:  1)  distance 

??«  si  "  DO,V'T-  measure  techniques  which  compare 

niny  min  infinite  resolution  at  low 

elevation  angle.  ground  truth  (or  model)  to  data  and 

minimize  a  suitable  measure  of  deviation 
between  the  two;  2)  correlation  techniques 
that  maximize  the  correlation  coefficient  between  the  model  pattern  and  data;  and 
3)  pattern  recognition  techniques  which  classify  data  according  to  a  set  of 
attributes. 

Distance  measure  techniques  require  a  suitable  definition  of  a  distance 
between  a  true  pattern  vector  and  the  data.  The  magnitude  of  this  distance  is 
then  compared  to  a  set  of  thresholds  which  are  boundaries  of  decision  regions. 

Two  distances  ha^e  been  considered: 

The  Euclidean  distance 


dR  --  jx-m|  (1) 

where  x  is  the  data  vector  and  m  the  true,  or  model,  vector  (Tou  and  Gonzales, 
1974).  For  instance,  m  could  consist  of  an  array  of  velocities  normalized  to  the 
maximum  (or  to  the  rms  value).  As  an  example,  nine  Doppler  velocities,  from 
three  consecutive  azimuths  and  three  consecutive  ranges,  of  an  ideal  mesocyclone 
centered  on  this  two-dimensional  array  could  form  the  mean  vector  rfi .  The  distances 


between  this  mean  vector  and  corresponding  data  vectors  must  be  compared  with  a 
suitable  threshold  and  when  they  are  smaller  than  this  threshold,  a  detection  of 
the  mesocyclone  is  declared. 

The  Mahalanobis  distance  (Tou  and  Gonzales,  1974) 


<L  =  (x-m)T  C'1  (141) 


(2) 


is  actually  the  negative  of  the  logarithm  of  the  likelihood  function  (if  x-m  has 


a  Gaussian  distribution).  C"^  is  the  inverse  of  the  covariance  matrix  of  the 


pattern  population  and  T  signifies  a  transpose  of  a  vector. 


Distance  measure  dm  requires  a  large  number  of  calculations  and  furthermore 


determination  of  mean  pattern  vector  m  and  the  covariance  matrix  C  is  complicated. 
The  complication  arises  because  the  mesocyclone  signature  of  a  vortex  with  a 
constant  diameter  is  a  function  of  range  and  at  a  same  range  mesocyclones  of 
different  size  would  require  mean  vectors  with  different  number  of  elements. 
Therefore,  an  adaptive  mean  vector  would  have  to  be  constructed. 

Correlation  techniques  have  pretty  much  the  same  type  of  problems  as  distance 
measuring  techniques.  They  require  maximization  of  a  correlation  coefficient 
between  the  mean  vector  and  the  data  vector. 

Substantial  storage  of  data  is  required  when  distance  or  correlation  techniques 
are  used.  If  the  pattern  vector  for  the  model  consists  of  data  from  M  radial, 
then  at  any  one  time  at  least  M  radials  of  data  must  be  stored  in  memory  so  that 
they  can  fast  be  made  available  for  calculations. 

The  pattern  recognition  method  employed  in  this  study  relies  on  recognition 
of  attributes  that  mesocyclones  possess.  Significant  attributes  have  been 
obtained  from  previous  measurements  (Burgess  et  al.,  1979).  Our  algorithm 
requires  simultaneous  storage  of  only  two  radials  of  velocities  and  should, 
therefore,  be  suitable  for  real-time  applications. 

Table  1 


Technique 


Number  of  Operations 


Number  of 
Storage  Locations 


Distance  df 


2M2-(N-M+1)2 


Distance  d, 


'm 


Correlation 
Pattern  Recognition 


2(2M2+M4)(N-M+1)2 

2M2(N-M+1)2 

M2(N-M+1)2 


M-N 

M-N 

M-N 

2N 


The  number  of  computations  and  storage  requirements  for  the  four  pattern 
recognition  techniques  are  compared  in  Table  1.  It  is  assumed  that  the 


mesocyclone  signature  (i.e.,  pattern  vector  m)  consists  of  MxM  model  velocity 

points  (at  M  discrete  azimuths  and  M  range  locations).  Let  the  sector  of  data 

that  encompass  the  storm  consist  of  an  NxN  array  (range,  azimuth)  of  velocities. 

The  number  of  operations  for  the  first  three  algorithms  is  the  total  number  of 

2 

additions  and  multiplications.  For  example,  one  dQ  calculation  requires  M 

6 

multiplications  (one  for  each  data  point  to  generate  normalized  velocities) 

2 

and  M  subtractions  (when  the  model  vector  is  subtracted  from  chosen  subset 

2 

of  data  values).  Since  there  are  (N-M+l)  distances  to  be  ca  lated,  the  total 

2  2 

number  of  additions  and  multiplications  is  2  M  (N-M+l)  .  It  hown  in  the  next 
section  that  for  our  algorithm  the  operations  are  comparisons  velocities 
between  two  adjacent  radial s  and  at  the  same  range  locations.  advantages  of 
this  algorithm  are:  1)  It  adapts  itself  automatically  to  the  size  of  the  mesocycl 
that  is,  there  is  no  need  to  assume  a  representative  size  as  ground  truth; 

2)  operations  are  comparisons  which  are  much  faster  than  multiplications  or 
additions;  3)  it  requires  storage  of  only  two  consecutive  radials  of  data  in 
addition  to  the  candidate  vectors  that  may  be  part  of  the  mesocyclone  or  shear 
lines;  and  4)  it  does  not  require  normalization  by  size  or  velocity. 

4.  THE  ALGORITHM 

The  first  step  in  the  algorithm  (flow  chart  in  Figure  2)  is  a  searcn  for 
floppier  velocity  gradients.  This  is  accomplished  by  identifying  either  a  con¬ 
sistent  increase  or  decrease  of  velocities  in  the  azimuthal  direction.  The 
algorithm  is  designed  to  recognize  only  cyclonic  rotation  and  was  meant  to  do  so 
cn  the  velocities  within  the  mesocyclone  core.  Therefore,  with  clockwise  antenna 
‘■otation  increasing  velocities  are  sought  and  with  counterclockwise  rotation  the 
decreasing  ones.  When  a  run  of  (at  least  two  data  points  are  required  to  form  a 
run)  increasing  or  decreasing  velocities  terminates,  a  pattern  vector  is  formed. 
The  five  components  of  the  pattern  vector  are  the  beginning  azimuth  <j>.  ,  the 
enemy  azimuth  tne  beginning  velocity  v^,  the  ending  velocity  v  ,  and  the 
'•ange  r.  Each  pattern  vector  must  pass  three  tests  before  it  is  saved  for  future 
analysis. 

The  first  test  the  vector  must  pass  is  a  low  "Doppler  angular  momentum" 
threshold.  The  pattern  vector  must  originate  from  data  that  have  sufficient 
"angular  momentum"  to  be  considered  for  further  analysis.  The  momentum  is 
estimated  from  the  product  of  the  measured  Doppler  velocity  difference  (v  -v.) 


START 
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and  the  azimuthal  distance  r(<J>  -<|>b)  stored  in  the  pattern  vector.  Therefore,  the 
calculated  momentum  often  may  not  correspond  to  a  rotating  air  mass,  but  rather 
to  a  shear  region.  Nevertheless,  this  first  step  is  designed  to  filter  out  a 
large  percentage  of  the  radar  data  so  that  the  more  extensive  analysis  is  performed 
only  on  the  remaining  pertinent  data.  A  mesocyclone  at  its  earliest  stages  may 
have  a  large  diameter,  and  while  the  shear  may  not  be  significant,  the  angular 
momentum  may  be  higher  than  it  will  be  throughout  the  rest  of  its  lifetime 
(i.e.,  the  angular  momentum  as  measured  by  the  radar  in  a  horizontal  plane  may 
decrease  towards  tornado  time).  The  second  test  uses  a  low  shear  threshold  to 
remove  anything  below  the  background  shear  value  of  2*10  s  (see  Burgess, 

1976).  In  the  third  test,  the  pattern  vector  is  reported  if  it  passes  either  a 
high  momentum  threshold  or  a  high  shear  threshold.  Since  shear  is  the  more 
dominant  feature  during  the  mature  stage,  it  is  u^ed  even  if  the  pattern  vector 
does  not  pass  the  high  momentum  threshold.  Two  thresholds  are  needed  to  accom¬ 
modate  variations  in  the  mesocyclone  from  its  early  stages  to  the  more  mature 
stage. 

Figure  3  shows  the  detection  region  in  the  momentum-shear  plane.  L  and  L 

s  m 

are  low  shear  and  momentum  thresholds  whereas  H  and  H  are  the  hi qh  thresholds. 

s  m  3 

Tests  on  data  are  needed  to  establish  optimum  values  of  these  thresholds,  and  if 
they  should  vary  with  range  from  radar.  In  the  present  study  the  following 
values  have  been  utilized: 

Ls=2  m-s  Vkm  Hs=4  m's~Vkm 
Lm=50  m-s”^ -km  H  =150  m-s”^ -km 

DETECTION  REGION 

Figure  Decision  boundary  in  the 
shear -women  turn  space .  The  thresholds 
_ ,  are  indicated . 
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An  alternate  view  of  the  detection  region  (Figure  4)  is  obtained  if  the 
velocity  difference  is  plotted  versus  the  azimuthal  distance.  In  both  Figure  3 
and  4  the  somewhat  peculiar  region  (l$  <  shear <  H$  and  Lm<  momentum  <  Hm)  of 
no  detection  is  clearly  visible.  This  whole  region  or  a  part  of  it  could  be 
included  into  the  detection  area  in  which  case  we  would  expect  an  increase  in 
both  the  probability  of  false  alarm  and  detection. 

After  a  field  of  data  has  been  filtered  by  the  aforementioned  tests,  the 
pattern  vectors  are  consolidated  according  to  their  relative  spatial  proximities. 
Azimuthal  and  range  distances  between  one  (comparing)  vector  and  all  other  already 
sorted  vectors  are  computed.  Whenever  the  comparing  vector  is  found  with  a 
center  closer  than  2.2°  in  azimuth  and  1  km  in  range  to  a  sorted  vector,  it  is 
put  in  the  same  feature  as  the  sorted  vector.  These  values  are  not  based  on  firm 
theoretical  or  experimental  grounds.  The  one  kilometer  separation  in  range  will 
insure  that  at  least  two  pattern  vectors  are  contained  within  a  rotation  that 
has  a  diameter  of  2  km  which  is  considered  to  be  about  the  minimum  for  mesocyclones. 
A  larger  spacing  in  azimuth  (for  ranges  over  30  km  from  the  radar)  allows  to 
group  vectors  that  are  azimuthal ly  overlapping  but  have  centers  that  are  further 
apart.  Optimum  value  of  the  azimuthal  separation  may  be  a  function  of  range  and 
together  with  the  range  separation  needs  to  be  further  investigated.  When  the 
comparisons  of  the  current  comparing  vector  with  sorted  vectors  is  exhausted,  one 
of  the  not  yet  sorted  vectors  becomes  a  comparing  vector  and  the  process  is 
repeated.  The  output  of  the  sorting  routine  is  a  set  of  features,  each  consist¬ 
ing  of  several  pattern  vectors.  A  feature  with  only  a  few  pattern  vectors  is 
discarded.  The  minimum  number  of  pattern  vectors  in  a  candidate  feature  must  be 
such  to  encompass  a  mesocyclone  with  minimum  core  diameter  (about  2  km).  The 
average  and  maximum  shears,  average  and  maximum  rotational  speeds,  average  momen¬ 
tum,  and  two  lengths  (diameters)  are  calculated  for  each  feature.  One  length 
is  for  the  radial  size,  and  the  other  for  the  azimuthal  size.  If  the  two 
lengths  disagree  by  more  than  a  factor  of  two  (i.e.,  the  feature  is  not  symme¬ 
trical),  the  feature  is  discarded.  However,  when  the  radial  length  is  large 
compared  to  the  azimuthal  one  or  vice  versa,  the  feature  can  be  classified  as  a 
shear  re g ion. 

iwo  slightly  different  methods  are  suggested  to  calculate  the  average 
azimuthal  "diameter".  In  one  the  "diameter"  D]  is  weighted  by  the  "momentum", 
and  is  calculated  from: 
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Figure  5.  Solidly  rotating  mesocyclonic  core.  Quantities  used 
in  the  text  to  define  various  calculations  of  diameters  are 
indicated. 

where  the  azimuthal  distance  =  r^.  at  an  i**1  range  gate,  v-j=vei~vbi ’  an(* 
i  bi  Tei 

For  a  combined  Rankine  model,  the  inner  core  of  the  mesocyclone  is  assumed 
to  be  a  solidly  rotating  cylinder  (Figure  5).  Then  the  Doppler  velocities  on  the 
periphery  of  the  cyclinder  are  given  by: 


v.,  =  -v_cosG;  v  .=vm  cosG  (4) 

bi  m  ei  m 

If  we  replace  the  summation  in  3  with  an  integral,  we  obtain 

tt/2  _  n/2  ? 

D-j  =  D  /  cosJ0/  /  COS^0=8D/3n  (5) 

O  0 

Therefore  in  order  to  obtain  the  true  diameter  D,  one  must  multiply  the  momentum 
weighted  diameter  with  3tt/8. 

The  second  method  calculates  the  "diameter"  from: 

D2  =  }>?/)>.  (6) 
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which  similarly  to  (5)  produces  a  theoretical  relationship: 
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■n/2  ?  rr/2 

D2  =  D  J  cos40d0/  /  cos0b0  *  Dtt/4  (7) 

o  0 

Factors  multiplying  0  in  (5)  and  (7)  are  close,  and  it  remains  to  be  determined 
which  weighting  defines  better  a  mesocyclone  diameter. 

A  correction  factor  would  also  be  needed  if  a  simple  average  of  azimuthal 
distances  is  used  for  estimating  the  diameter.  The  relationships  then  are: 

*/2 

03  =  =  2  /  Dcos0=2D  (8) 

0 

Our  experience  with  mesocyclones  is  that  the  Doppler  velocity  couplet  defining 
the  maximum  speeds  of  mesocyclones  (smallest  closed  contours  on  Figure  1)  conforms 
well  with  theoretical  predictions  but  other  velocity  contours  may  deviate  consider¬ 
ably  from  the  model.  Weighting  by  the  momentum  restores  some  symmetry  and  gives 
significance  to  the  data  that  inherently  correspond  better  to  the  model.  Further¬ 
more  data  across  the  maximum  Ai  would  be  less  contaminated  by  errors.  For  these 
same  reasons,  azimuthal  center  <{>  and  range  center  rc  of  the  mesocyclone  are 
obtained  by  weighting  the  ranges  r^  and  azimuthal  centers  <t>  .  =  (<P^j+<Pe^  )/2  of 
vectors  wi th  ^ 


rc =  IViVlVj  (9) 

*c  =  (,0) 


Calculations  of  and  in  the  program  are  slightly  at  variance  with  (3)  and 
(6).  In  order  to  preserve  computational  speed,  r^  which  is  almost  constant  is 
nulled  out  of  the  summation  and  replaced  with  (9),  thus: 


D1  =  rc  friV^Vi 
D2  =  rc  ^i  /  ^i 


(11) 

(12) 
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NORMALIZED  DISTANCE  BETWEEN  BEAM  CENTER  AND  VORTEX  CENTER  x0/r, 


Figure  6.  Normalized  Doppler  velocities  of  a  combined  Rankine  model  vortex 
centered  at  a  radar  range  r.  Vortex  radius  is  rt  and  the  range  resolution  is  much 
finer  than  the  size  of  the  mesocyclone . 


Presently,  the  program  classifies  as  mesocyclones  all  features  for  which  the 
ratio  of  the  radial  to  azimuthal  diameters  lies  within  0.5  to  2.  We  need  to  gain 
more  experience  in  order  to  firm  up  this  criterion.  However,  consideration  of 
resolution  at  the  furthest  range  of  interest  (250  km)  points  out  that,  with  a  1° 
antenna  beamwidth,  a  mesocyclone  having  a  solid  core  diameter  of  2  km  (out  of  37 
mesocylones  observed  between  1971  and  1975,  only  in  one  case  was  the  diameter  as 
small  as  2  km  (Burgess,  1976))  would  seem  to  have  an  azimuthal  diameter  twice  as 
large  (i .e. ,  4  km) . 

This  can  be  demonstrated  observing  Figure  6  which  is  a  plot  of  normalized 
mean  Doppler  velocity  versus  normalized  distance  between  beam  center  and  vortex 
center  for  a  radar  range  that  is  centered  on  the  vortex.  The  vortex  is  a  Rankine 
combined  model  and  range  resolution  is  assumed  to  be  perfect. 

5.  PRELIMINARY  RESULTS 

To  date,  we  have  preliminary  results  from  two  cases.  The  first  storm 
produced  the  mesocyclone  on  April  30,  1978,  and  spawned  a  maxi-tornado  which 
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destroyed  25  homes  in  and  near  Piedmont,  Oklahoma  (Burgess  et  al .  1979).  Five 
different  times  have  been  examined. 

Comparisons  between  mesocyclone  parameters  obtained  by  the  algorithm  and 
subjectively  deduced  values  from  mean  Doppler  velocities  are  shown  in  Figure  7 
and  Figure  8.  The  agreement,  in  general,  is  quite  good;  for  all  five  locations 
the  two  positions  are  within  the  radius  of  the  mesocyclone  and  the  tracks  are 
almost  identical.  Note  that  the  subjectively  determined  diameters  are  somewhat 
smaller.  This  is  not  surprising  because  the  two  methods  for  estimating  the 
diameter  differ.  The  meteorologist  who  estimated  the  diameters,  used  the  distance 
between  peak  velocities  (i.e.,  the  farthest  closed  isodops  as  in  Figure  1).  On 
the  other  hand,  the  algorithm  calculates  the  weighted  average  of  distances  between 
the  peak  velocities  in  a  feature  (Eq.  11).  No  attempt  has  yet  been  made  to 
account  for  smoothing  by  the  radar  resolution  volume. 

The  average  shear  calculated  by  the  algorithm  from  the  pattern  vectors  of 
the  mesocyclone  feature  follows  the  same  trend  as  the  maximum  shear  determined 
subjectively  (Figure  8).  A  somewhat  large  difference  occurs  at  the  latest  time 
when  the  maxi -tornado  formed  and  the  subjectively  measured  shear  was  influenced 
by  the  tornado  vortex  signature.  Since  the  maximum  shear  is  contained  in  one  of 
the  pattern  vectors,  with  a  simple  modification  of  the  program,  it  could  be 
easily  retrieved. 

The  second  storm  on  which  the  algorithm  was  tested  occurred  on  May  22,  1981, 
west  of  Norman  and  produced  several  very  large  tornadoes.  One  example  of  a 
mesocyclone  feature  in  a  range  versus  azimuth  display  of  velocities  is  shown  on 
Figure  9.  The  total  of  32  vectors  are  in  this  feature.  Two  vectors,  bracketed  by 
dasned  lines  towards  the  top, were  not  detected  because  of  weak  momentum. 

However,  the  vector  above  them  was  detected  and  placed  in  the  feature  due  to  its 
proximity  to  the  vector  at  68.68  km.  The  star  on  this  figure  locates  the  meso¬ 
cyclone  center  computed  from  Eqs.  (9),  (10). 

Figure  10  illustrates  what  appears  to  be  two  separate  features  detected  by 
the  algorithm.  Upon  closer  examination  of  the  velocity  field,  we  note  that  there 
were  several  non-valid  velocities  (999)  between  the  two  features.  The  algorithm 
would  have  detected  four  vectors  (at  65.8,  65.9,  66.1,  66.2  km)  if  the  non-valid 
velocities  at  282.7°  and  281.2°  were  bypassed.  Then  the  gap  between  the  closest 
vectors  in  range  would  have  been  less  than  1  km  and  all  would  be  classified  as  a 
single  feature.  The  following  modification  could  accommodate  such  situations: 

If  there  is  an  isolated  (in  azimuth)  non-valid  data  (i.e.,  velocity  at  a  very  weak 
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Fixture  9 .  Mesocuclon e  feature  detected  by  the  algorithm .  Velocities  are  printed 
in  an  azimuth  (on  abscissa)  versus  range  7 ordinate )  format .  Elevation  angle 
is  7.*?°.  Solul  Iznes  bracket  all  detected  vectors  of  a  feature .  Dashed  lines 
Carres  rend  to  two  nondetected  vectors . 


reflectivity  or  a  radial  of  bad  values),  assign  to  it  the  value  of  the  adjacent 
(previous  azimuth)  velocity  and  proceed  with  normal  vector  calculations.  If 
there  are  2  bad  adjacent  velocities,  stop  the  run  and  start  the  new  vector.  In 
this  way,  still  only  two  radials  at  a  time  are  needed  and  one  additional  array 
of  flag  values,  and  yet  isolated  bad  data  would  not  cause  disruption  of  the 
vector  calculations.  Of  course,  at  the  expense  of  more  storage  and  computing 
time  one  could  treat  two  consecutive  bad  data  as  well  as  situations  where  there 
is  a  slight  reversal  in  gradient  (i.e.,  one  or  two  data  points).  It  is  note¬ 
worthy  that  the  vector  at  70.03  km  ( v^=- 51  m-s"\  vg=26  m-s"^)  was  not  detected. 
Even  though  the  shear  is  very  large,  the  momentum  is  weak  and  thus  the  data  do 
not  pass  our  test.  Because  such  situations  may  occur  with  tornadoes,  it  may 
prove  useful  to  include  an  extremely  high  shear  threshold  in  combination  with 
our  test. 

Three  consecutive  scans  (in  elevation)  with  the  size  and  location  of  a 
mesocyclone  detected  by  the  algorithm  are  shown  on  Figure  11a,  b,  and  c.  Note 
again  a  split  detection  of  the  mesocyclone  at  times  1  and  2  (Figure  11a).  Tne 
case  at  time  1  corresponds  to  Figure  10,  but  the  split  at  time  2  occurred  for  a 
different  reason:  Namely,  the  momentum  and  shear  between  the  two  features  were 
very  weak. 

On  the  scan  above  (Figure  lib)  we  note  another  split  signature  at  time  2. 

The  reason  for  the  split  now  is  that  the  vector  centers  are  separated  in  azimuth 
(although  close  in  range).  It  appears  that  intense  circulation  is  developing  at 
65  km  west  and  9  km  north  (A)  where  the  maximum  shear  was  55  m-s"Vkm  and  maximum 
momentum  is  405  m-s'^-kni.  The  feature  whose  center  is  marked  with  o  had  the 
following  maximum  shear  and  momentums:  at  time  1,  max  shear  =  9  m-s"Vkm; 
max  momentum  =  405  m-s'^-km.  Clearly  the  A  feature  is  part  of  the  larger  circu¬ 
lation  and  was  contained  in  its  feature  at  time  1  as  evidenced  by  the  large 
'hear.  With  the  split  the  feature  A  took  with  it  the  large  shear.  The  split 
here  could  have  been  prevented  if  somewhat  larger  distances  between  vectors  were 
allowed  in  the  sorting  routine.  This  criterion  is  again  a  subject  of  current 
investigations . 

The  Figure  11c  at  5.4”  elevation  shows  two  consecutive  positions  of  the 
mesocyclone  (time  1  and  2)  which  were  quite  well  detected. 


16 


NORTH  (km)  NORTH  (km) 


10 


Figure  11a .  The  position  and  sizes  of 
features  detected  at  two  consecutive 
scans  and  same  elevation  3.3°.  (1  and 
2)  separated  by  6  minutes .  Azimuthal 
diameters  are  indicated  with  vertical 
lines  and  range  diameters  with  hori¬ 
zontal  lines .  Two  features  A  and  0 
were  detected  in  both  scans .  The  two 
features  in  the  first  scan  (1)  cor¬ 
respond  to  the  ones  on  Figure  9. 


Figure  lib .  Same  as  11a  but  at  the 
next  elevation  angle  (4.5°)* 
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6.  CONCLUSION 

A  novel  pattern  recognition  algorithm  for  detecting  mesocyclonic  shear 
from  single  Doppler  radar  data  has  been  developed.  The  method  is  based  on  extract¬ 
ing  features  that  characterize  mesocyclones.  It  operates  simultaneously  on  two 
radials  of  velocities,  utilizes  relatively  simple  computational  techniques,  and 
should,  therefore,  be  applicable  to  real-time  operational  systems.  Preliminary 
resu  ts  on  a  limited  data  set  are  encouraging;  however,  further  comprehensive 
test  .'i..  is  needed  to  establish  optimum  thresholds  and  their  range  dependence. 

Dd'*ing  the  course  of  this  study,  some  simple  modifications  of  the  algorithm 
rhit  could  increase  its  accuracy  and  usefulness  became  apparent.  1)  For  instance, 
the  maximum  measured  shear  and  momentum  could  be  printed  out.  2)  In  order  to 
ivoid  loss  OT  vectors  with  very  high  shear  but  small  momentum  such  as  in  the 
-enter  of  the  tornado  vortex  signature,  it  may  be  advisable  to  include  in  an  OR 
combination  with  our  test  a  very  high  shear  test.  3)  Presently,  an  isolated  bad 
velocity  (in  azimuth)  signals  the  end  of  a  vector;  this  need  not  be  so,  to  pre¬ 
serve  storage  and  speed  two  consecutive  invalid  points  could  signal  the  end  of  a 
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vector.  4)  The  optimum  values  of  azimuthal  and  range  distance  that  vectors  must 
satisfy  to  form  a  feature  should  be  established.  5)  Optimum  momentum  and  shear 
thresholds  need  to  be  determined. 

The  present  version  of  the  algorithm  does  not  consider  the  time  and  height 
continuity  of  the  circulation  but  such  an  addition  is  quite  straightforward. 

There  may  be  other  modifications  that  we  cannot  foresee  without  extensive  tests 
on  various  data  sets.  This  would  also  help  us  establish  the  false  alarm  rates  and 
detection  probabilities  per  scan. 
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APPENDIX  A 


Program  Description 

The  mesocyclone  detection  program  requires  radials  of  preprocessed 
velocity  data.  Preprocessing  involves  assignment  of  correct  ranges  to 
velocity  estimates,  thresholding  and  properly  scaling  these  estimates,  and 
reading  and  sorting  out  various  housekeeping  information.  Preprocessing 
the  mesocyclone  program  and  the  various  subroutines  of  the  mesocyclone  program 
are  described  in  this  appendix. 

A1 .  Preprocessing 

Preprocessing  starts  with  the  routine  RDREC  that  is  called  for  every 
radial  (see  Figure  A.l). 

For  Every  Radial : 

1.  Read  raw  data  and  strip  and  identify  housekeeping--subroutine  RWDP81 
Sample  call:  CALL  RWDP81 (W,  ISKP,  LTIM,  LINT) 

W--logical  unit  number  assigned  to  tape  drive 
ISKP-1:  skip  1  record 

2:  decode  housekeeping  only 

3:  decode  housekeeping,  intensities  and  time  data 
4:  decode  intensities  and  time  series  without  reading  a  record 
800:  change  mode 

0:  reset  all  gate  alignments  to  1 
9:  rewind  tape 

LTIM- 1 :  decode  time  series  data 
0:  don't  decode  time  series 
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LINT-1:  decode  integrator  PPP  data 
0:  don't  decode  PPP  data 

The  typical  call  from  the  mesocyclone  program  is:  CALL  RWDP81 (10,3,0,1 ) 

Assign  correct  ranges  to  velocity  estimates.  This  is  needed  when 
interlaced  transmission  for  estimating  separately  velocities  and 
reflectivities  is  used. 

Subroutine  DPEXPAND 

Sample  call:  CALL  DPEXPAND(IVCINT,RVCPPP,RSDPPP,STN,IB,IE,IDEL, 
SNRT,VT,STDT,INT,VEL,STD) 

IVCINT-integrator  array  from  RWDP81 
RVCPPP-velocity  array  from  RWDP81 
RSDPPP-spectrum  width  array  from  RWDP81 

STN-signal-to-noise  ratio  array  -  64  values 
IB-beginning  range  gate  for  processing 
IE-ending  range  gate  for  processing 
IUEL-requi red  gate  spacing 
SNRT-si gnal -to-noise  ratio  threshold 
7T-velocity  threshold 
STDT-spectrum  width  threshold 
INT-expanded  integrator  array 
VEL-deal iased  velocity  array 

STD-unbiased  spectrum  width  array  (not  used  by  the  program) 

Typical  call  from  mesocyclone  program:  CALL  DPEXPAND ( IVCINT,  RVCPPP, 
RSDPPP,  STH, 1,1524, 1,0, 10, 15, OUT,  0UT1  ,  0UT2) 

where  OUT  is  an  integar  array  dimensioned  1524  and  0UT1  and  0UT2  are 
real  arrays  dimensioned  1524 


At  NSSL  the  interlaced  transmission  uses  a  pulse  separation  time  for 
reflectivity  estimation  that  is  four  times  the  one  for  velocity  estimation. 
However,  in  our  program  only  the  first  two  ambiguous  trips  for  velocity 
estimation  are  considered.  Since  the  number  of  usable  range  gates  for  velocity 
estimation  is  762, for  two  trips  the  total  number  of  locations  is  twice  762 
or  1524. 

3.  Calculate  reflectivities:  subroutine  CALREF 

Sample  call:  CALL  CALREF  (OUT,  REFL,STN,  RANGE,  LOGR, CNST) 

OUT:  integrator  array  passed  from  DPEXPANO 
REFL:  array  of  reflectivities  at  each  range  gate  processed  in 
DPEXPAND  (dimensioned  same  as  OUT) 

STN:  signal-to-noise  ratio  array  (64  values) 

RANGE:  array  (dimensioned  same  as  OUT)  containing  the  range  in  kilometers 
to  each  range  gate 

LOGR:  array  (dimensioned  same  as  RANGE)  containing  L0GR( I )=2Log^ 
(RANGE(I))  where  I  represents  the  range  gate  number 
CNST:  a  constant  needed  for  calculation  of  reflectivities 

CNST=CN*1 0- 110.8  where  CN  is  the  result  of  a  call  to  the 
function  RADARCON 

RANGE,  LOGR,  and  CNST  are  calculated  in  subroutine  PRELIM  of  the 
mesocyclone  program. 

The  call  of  subroutine  CALREF  in  the  program  looks  exactly  like  the 
sample  call  above. 

For  E a ch  Sc an : 

1.  calculate  range  to  the  first  range  gate  and  distance  between  range 
gates  in  kilometers  -  subroutine  D0PRNG81 
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Sample  call:  CALL  D0PRN68 1 ( DAT E , BR6 , DRG , V ) 

DATE:  day  data  were  collected  (i.e.,  052281) 

BRG:  range  to  first  range  gate  in  km 
DRG:  range  between  gates  in  km 
V:  Nyquist  velocity  (not  used  in  program) 

2.  Calculate  a  constant  (CN)  which  determines  the  constant  (CNST),  used  to 
calculate  reflectivities  -  function  RADARCON 
Sample  call:  CN= RADARCON (BH,BW,P,WL,G,P0W) 

BH:  vertical  beamwidth  in  degrees 
BW:  horizontal  beamwidth  in  degrees 
P:  pulse  width  in  microseconds 

WL:  wavelength  of  transmitted  electromagnetic  field  in  centimeters 
G:  antenna  gain  in  DB 
POW:  peak  transmitter  power  in  kilowatts 
Typical  call  from  mesocycl one  program:  CN=RADARC0N(0. 8, 0.8,1 .0,10.52, 
46.8,750.0) 
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Mesocyclone  Program 


Input  to  Mesocyclone  Proaram 

From  preprocessing  -  passed  through  common  blocks  and  subroutine  calls 

1.  one  radial  of  data 

OUT!  -  array  returned  from  subroutine  DPEXPAND  contains  real 

velocities,  which  are  changed  to  integers  and  assigned  to 
array  VEL  by  subroutine  RDREC  in  meso  program. 

REFL  -  array  returned  from  subroutine  CALREF  -  contains 

reflectivities.  These  are  checked  against  a  lowest  needed 
reflectivity  value  (RTH),  and  velocities  corresponding  to 
reflectivities  below  this  level  are  not  processed. 

2.  Housekeeping  Data 

NTIME  -  time  current  radial  was  collected 
AZM  -  azimuth  angle  of  current  radial 
ELEV  -  antenna  elevation  angle  during  current  scan 
STATUS  -  tape  status  indicator  (after  last  read) 

0:  end-of-file 
1 :  bad  record  was  read 
-1:  good  record  was  read 
ITYPE  -  record  type 

0:  Bad  1:  Time  series  2:  Pulse  pair 
ISTAT  -  Station  code 


1:  Norman  2:  Cimarron 


IPRTSW  -  pulse  repetition  time  indicator 


RPRT  -  constant  proportional  to  the  pulse  repetition  time. 
The  relationship  between  the  switch  and  RPRT  and  the 
true  pulse  repetition  time  is  as  follows: 


IPRTSW 

RPRT 

True  pulse  repetition  time 

RPRT*0. 768  ms 

0 

1.0 

0.768 

ms 

1 

1.2 

0.921 

ms 

2 

1.4 

1.075 

ms 

3 

1.6 

1.228 

ms 

IDELAZ  -  azimuth  step  value  indicator 
1:  0.50  2:  1.0°  3:  2.0° 


IAZMD  -  azimuthal  antenna  direction  indicator 
1:  counterclockwise  3:  clockwise 
B.  OTHER  INPUT  DATA 
1.  read  from  files 

SHEART,HSHEAR, THRESH, LOWTH  -  thresholds  used  in  detecting  velocity 
runs.  See  explanation  of  subroutine  SELECT 

M,N,R  -  values  used  to  sort  features  -  see  description  of  subroutine 
FEAT 

RTH  -  reflectivity  threshold  -  it  has  been  found  experimentally 

that  mesocyclones  do  not  occur  in  regions  with  reflectivities  below 

this  threshold.  Therefore,  to  save  computer  time,  velocities 

corresponding  to  reflectivities  below  this  level  are  ignored. 

ELTH  -  elevation  limit  -  This  limit  is  set  to  correspond  to  a 

height  of  about  8  km  at  the  storm  location.  Above  this  height  the 

circulation  signature  changes  into  a  divergence  pattern  and  can  not 

be  detected  with  the  mesocyclone  algorithm. 
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OUTPUT  DATA 

One  line  of  output  is  printed  at  the  beginning  of  each  scan.  Also 
each  feature's  vectors  are  printed  out  along  with  feature  attributes  after 
each  scan. 

A.  Printouts  at  the  beginning  of  scan  are: 

1.  Direction  of  the  antenna  in  coded  form 
(see  IAZMD  under  housekeeping  input) 

2.  Elevation  angle  of  the  radar  antenna  (see  ELEV  under  housekeeping 
input) 

3.  The  time  of  the  collection  of  the  first  radial  of  the  scan 
(see  NTIME  under  housekeeping  input) 

OUTPUT  FOR  EACH  FEATURE: 

B.  Features  -  array  L I  ST ( I  ,J,K)  where  I  =  feature  number,  J  =  vector 
number  and  K  =  attribute  of  vector  number. 

A  feature  is  a  group  of  vectors,  each  vector  representing  a  velocity 
run  at  a  constant  range.  Each  vector  has  the  following  seven  attributes 

1.  range  to  the  run  in  kilometers 

2.  azimuth  angle  where  the  run  began 

3.  azimuth  angle  where  the  run  finished 

4.  velocity  recorded  where  the  run  began 

5.  velocity  recorded  where  the  run  ended 

6.  the  calculated  momentum  of  the  run. 

7.  the  calculated  shear  of  the  run. 


Feature  attributes  -  array  AFEAT(I.J)  where  I  =  feature  number, 

J  =  attribute  (of  feature)  number 

1.  number  of  vectors  in  the  feature 

2.  lowest  azimuth  angle  in  the  feature  in  degrees 

3.  highest  azimuth  angle  in  the  feature  in  degrees 

4.  lowest  range  in  feature  in  kilometers 

5.  highest  range  in  feature  in  kilometers 
Other  feature  attributes  (not  in  the  array  AFEAT) 

1.  momentum  weighted  center  of  the  feature  (azimuth  and  range,  array  PT) 

2.  azimuthal  diameter  of  feature  -  weighted  by  momentum  (OIAMAZ) 

3.  unweighted  azimuthal  diameter  (ODIAM) 

4.  unweighted  radial  diameter  (DIAMR) 

5.  ratio  DIAMAZ/DIAMR  for  symmetry  considerations 

6.  Average  shear  of  feature  (AV3SHR) 

7.  Maximum  shear  of  feature  (MAXSH) 

8.  Average  momentum  of  the  feature  (AVGMOM) 

9.  Average  rotational  speed  of  feature  (AVGSP) 

0.  Maximum  rotational  speed  of  feature  (MAXSP) 


A. 3.  Description  of  Routines 
MAIN  PROGRAM 

The  main  program's  function  is  to  control  the  processing  done  by  the 
subroutines.  Its  flow  chart  is  on  Figure  A. 1  and  the  control  proceeds  in 
the  following  order: 

Prel i mi  nary 

1.  read  variables  from  files 

2.  initialize  STATUS, NTIME, DONE 

3.  enter  beginning  and  ending  processing  times  on  terminal 
( IBEG, IENDT) 

For  Every  Scan 

4.  zero  out  vector  and  feature  lists  (VECTOR, LIST) 

5.  initialize  ERR0R(=0),C0UNT(=0),  and  FIRST(=1) 

6.  read  one  needed  radial  (CALL  RDREC) 

7.  check  STATUS  0:  end  program  else  continue 

8.  write  important  information 

9.  zero  out  CODE  array  and  reassign  VEL  to  OLDVEL 

10.  dealias  OLDVEL  and  reassign  housekeeping 

11.  check  DONE  1:  end  program  else  continue 

12.  read  one  needed  radial  (CALL  RDREC)  (FIRST=0) 

13.  check  STATUS  0:  D0NE=1 

14.  bad  radial  yes:  Go  to  11  No:  continue 

15.  end  of  scan?  yes:  reassign  aximuthal  direction  (IAD)  and  go  to  19 
no:  continue 

16.  dealias  VEL,  search  for  runs  (CALL  SELECT) 

17.  0LDVEL=VEL,  reassign  housekeeping 


18.  go  to  11 

19.  VEL=999's  (means  not  usable  velocities),  reassign  IAD 

20.  search  for  runs 

21.  no  vectors?  Yes:  go  to  23  No:  continue 

22.  Assemble  features  (CALL  FEAT)  and  calculate  feature  attributes 
(CALL  CALC) 

23.  check  DONE  0:  go  to  4  1:  end  program 


SUBROUTINE  SELECT 

The  flow  chart  of  this  subroutine  is  shown  on  Figure  A.2.  It  operates 
on  two  consecutive  radial s  of  velocity  data,  OLDVEL  and  VEL.  The  purpose 
of  the  routine  is  to  search  at  every  range  gate  for  changing  velocities  with 
respect  to  azimuth  angle.  The  array  CODE  keeps  track  of  velocity  trends 
(at  every  range  gate)  as  of  the  last  call  to  SELECT 

Values  of  CODE  elements  and  their  meanings 
Value  Tendency 

0  constant  velocities  (no  run  is  going  on) 

-1  decreasing  velocities  (run  is  on) 

1  increasing  velocity  (run  is  on) 

Ire  routine  executes  a  DO  loop  with  an  index  from  1  to  1524  called  I.  On 

each  execution  of  the  statements  in  the  loop,  OLDVEL(I)  is  compared  to  VEL(I). 

This  tendency  is  compared  to  CODE ( I ) .  999's  cause  a  termination  of  an  ongoing  run 
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Fiyuv<j  A.?.  Flow  chart  of  the  Subroutine  SF.bF.CT 
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0 


start  increasing  run 


constant  velocity 


start  decreasing  run 


1  continue  increasing  continue  increasing  finish  increasing 

run  run  run  and  start  decreasing 

run 

-1  continue  decreasing  finish  decreasing  run  continue  decreasing  run 

run  and  start  increasing 

run 

Control  proceeds  as  follows: 

1.  1  =  1 

J 2 .  a)  constant  velocities  or  continuing  run:  go  to  3 

b)  starting  a  run:  go  to  START  NEW  RUN 

c)  ending  a  run:  go  to  END  RUN 

3.  1=1+1 

4.  I  ’524:  return  to  main  program  else:  go  to  2 

END  Ri'N: 

a)  ENDV EL- { OLDV EL ( I )*1 .07)1 RPRT  last  velocity  of  run 
ENDAZ=01.DAZ  last  azimuth  of  run 
D£LVEL=ENDVEL  -  BEGBEL ( I )  velocity  change  of  run 
DELAZ- j ENDA7  -  BEGAZ(I)j  azimuth  change  of  run 

b)  Throw  out  unwanted  runs  (only  runs  consistent  with  cyclonic 
circulation  are  kept.) 

clockwise  antenna  rotation:  throw  out  if  DELVELiO 
counterclockwise  rotation:  throw  out  if  DELVELiO 
if  run  is  thrown  out:  go  to  START  NEW  RUN 
else:  continue 
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c)  LENGTH=DELAZ  CON  RANGE ( I )  C0N="/180 

DELVEL= | DELVEL | 

if  LENGTH  i  0:  go  to  START  NEW  RUN;  else:  continue 

d)  AMDM-DELVEL (LENGTH )  momentum  of  run 

SHR=DELVEL/LENGTH  shear  of  run 

e)  if  AMOM  >  LDWTH  and  SHR  >  HSHEAR  then  go  to  f) 

if  AMOM  ^  THRESH  or  SHR  <  SHEART  then  go  to  START  NEW  RUN 
else:  continue 

f)  save  vector  (C0UNT=C0UNT+1 ) 

too  many  vectors:  return  to  main  program 
else:  go  to  START  NEW  RUN 
START  NEW  RUN: 

a)  V EL ( I )=999 :  CODE ( I )-0,  go  to  e) 
else:  continue 

b)  increasing  velocities:  C0DE(I)=1,  go  to  e) 
else:  continue 

l)  decreasing  velocities:  CODE ( I )  =  - 1 ,  go  to  e) 
else:  continue 

d)  C0DE( I  )=0 

e)  BEGV EL ( I ) -OLDVEL ( I ) ( 1 .07)/RRl 
BEGAZ( I )=0LDAZ 


go  to  3. 
THRESHOLDS : 


Shear  and  momentum 
'  if.  shear: 


"•  . ’ i  shear: 


'  aw  momentum: 

,  ■  s";  momentum: 


thresholds  and  their  presently  used  values 
SHEART  =  2.0  (m- s' 1 'km) 

HSHEAR  =  T.O  (m- s' 1 /kin) 

LOWTH  ■  bO.O  (in -s'1  -km) 

THRESH  '  >0.0  (m-s'^k:-' 


(i  ro 
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SUBROUTINE  FEAT 


This  routine  (figure  A. 3)  sorts  vectors  found  by  SELECT  into  features 
which  could  be  mesoc.yclones ,  shear  lines  or  have  small  size  and  therefore  be 
rejected.  Every  vector  found  by  SELECT  during  one  scan  of  the  radar  is  assigned 
to  a  feature,  so  many  small  features  will  be  thrown  out  by  subroutine  CALC. 

Control  proceeds  as  follows: 

1.  Put  1st  vector  into  1st  feature.  Initialize  1st  feature's  attributes. 
Mark  (flag)  1st  vector  used. 

2.  COOL  =  0 

3.  If  no  more  vectors  -  go  to  7  else:  get  next  vector  (we'll  call  it  VA) 

A.  If  VA  is  marked  used,  go  to  3  else:  continue 

b.  If  no  more  vectors  in  current  feature,  then  go  to  3  ;  else  get 

current  feature  element  in  current  feature  (VB) 

6,  It  the  azimuthal  center  of  VA  is  within  2.2^  of  the  azimuthal  center 
of  VB  AND  the  range  of  the  2  vectors  are  within  1  km 

a)  put  VA  in  VB  *  s  feature 

b)  adjust  feature  attributes 

c)  mark  VA  used 

d)  COUL-l 

e )  Go  to  3 
else  go  to  b 

7.  It  C0:.l  1,  cm  to  2  (on  step  3  VA  =  first  vector)  else:  continue 
Got  -if’.;  unmarked  vector  (VA).  if  all  vectors  used,  go  to  1! 

:  1  <  '  r  ,<  new  feature  (INIJEX-INDEX+1 );  too  many  features,  go  to  11 
'r;ititl..e  feature  attributes 
G,  to  2  (on  ■  VAr first  vector) 

:■!.  ; .V;L  X  I NULX-  1 
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SUBROUTINE  CALC 

This  routine  finds  important  feature  data  used  to  discern  which  features 
rt,t*  lliesocyclones,  which  are  shear  lines,  and  which  are  neither.  The  control 
pattern  below  is  repeated  for  each  feature  found  in  the  last  scan,  using  a  do 

loop  with  index  !  (‘igure  A. 41. 

1.  Not  enough  vectors  in  feature  -  Go  to  11  else:  continue 

2.  NSIIM  =  «  vectors  in  feature.  Zero  out  all  sums 

J.  Compute  sums  using  do  loop  with  index  J.  The  sums  range  from  1  to 
NSUM. 

a}  Vector  attributes  used  (vectors  are  part  of  features) 

»  -  range:  LIST(1,J,1) 

BAZ  -  beginning  azimuth:  LIST(I,J,2) 

E.AZ  -  ending  azimuth:  LIST([,J,3) 

6V  -  beginning  velocity:  LIST(I,J,4) 

EV  -  ending  velocity:  LIST(i,J,5. 

MG  -  "momentum" ;  here  is  a  product  of  azimuthal  separation  (,;) 

with  velocity  difference:  LIST(I,J,6) 

SHEAR  -  shear=  |cV-3V|*l80/(DcLAZ*R*Jl) : LIST ( I , J , 7 ) 

MAZ  -  center  azimuth  *■  (EAZ-BAZ)/2+BAZ 
DLLAZ  -  change  in  azimuth  =  EAZ  -  BAZ 

Provision  is  ruade  to  account  for  the  cases  where  the  zero  azimuth 
is  between  the  beginning  and  ending  azimuths. 

b)  SUMS 


SUM-  )MC' 

SUMR  yp-MC 

SUMAZ7MAZ-m 


azm=j;deiaz-m0 

oaz=£delaz 

4.  Compute  data 

a)  AVGSHR  -  average  shear  of  feature  =  SHSUM/NSUM 

b)  AVGMOM  =  average  momentum  =  SUM/NSUM 

c)  DIAMR  -  radial  diameter  =  AFEAT( I ,5)-AFEAT( I ,4) 

d)  MAXSH  -  maximum  shear 

e)  MAXSP  -  maximum  rotational  speed  =  | EV-BV | /2 

5.  If  GUM  =  0  or  DIAMR  =  0  or  OAZ  =  0 
go  to  11 

else  continue 

6.  Compute  more  data 

a)  PT  -  center  of  feature  PT (1,1)  radial  =  SUMR/SUM 

PT(I,2)  azimuthal  =  SUMAZ/SUM 

7.  If  PT (1,1)  =  0 
go  to  11 

else  continue 

8.  Compute  more  data 

a)  QIAMAZ  -  azimuthal  diameter  weighted  by  momentum 

=  AZM- PT(  1 , 1 )  -  ii2- 3/  (8- 1 80- SUM) 

b)  ODIAM  -  azimuthal  diameter  -  unweighted 

-  AZM- PT (1,1 )/ (OAZ- 45) 

c)  RATIO  =  DIAMAZ/DIAMR  ) 

)  to  determine  symmetry 

d)  R2  -  0DTAM/DIAMR  ) 

e)  AVGSP  -  average  rotational  speed  =  AVGSHR*DIAMAZ/2 

9.  If  feature  is  symmetric,  it  is  a  mesocyclone  (.5- RATIO-,?  or  ,5>-R?  2) 
else  feature  is  a  shear  line 
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Flow  -hart  of  the  Subroutine  CALC. 
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10.  Write  feature  and  its  attributes 

11.  Stop. 


SUBROUTINE  RDREC 

RDREC  gets  records  from  tape  in  a  quasi -real  time  manner--one  record 
at  a  time.  Also  the  routine  calls  other  needed  preprocessing  routines  (Figure  A. 5) 
Control  proceeds  as  follows: 

1.  If  end-of-file  or  end-of-processing,  go  to  8. 

Else  continue 

2.  Read  1  radial  of  raw  data  (CALL  RWDP81 ) 

3.  Undesired  radial? 

STASUS  =  1 ,  bad  radial 

ITYPE  t  2,  undesired  type  of  data 

IAZMD  =  2,  stationary  antenna  (azimuthally) 

ISTAT  f  1,  Cimarron  data 

NTIME  <  IBEGT,  before  desired  processing  times 
ELEV  ?  ELTH,  too  high  for  mesocyclones 
Go  to  1,  else  continue 

4.  If  first  =  1,  call  PRELIM  and  continue,  else  continue 

5.  Expand  data  -  CALL  DPEXPAND 

6.  Calculate  reflectivities:  CALL  CALREF 

7.  If  reflectivity  at  any  gate  reflectivity  threshold,  then  ignore 
velocity  at  that  gate.  Go  to  9. 

8.  STATUS  -  0  Set  end-of-file  indicator 

9.  Return 
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RDREC 

COND  I  -  STATUS  *  0  AND 
NTIME  <  IENDT 

COND  2  *  STATUS  *  I  OR 

1TVPE  #  2  OR  IAZMD‘2 
OR  ISTAT  *  I  OR 
NTIME  <  1BEGT  OR 
ELEV>ELTH 
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SUBROUTINE  PRELIM 

PRELIM  calculates  the  pulse  repetition  time  and  the  range  to  each  range 
gate.  ALso  the  routine  finds  the  1 og^ Q  of  the  square  of  the  range  at  each 
range  gate  (Figure  A. 6).  Control  proceeds  as  follows: 

1.  Calculate  BRG  and  DRG  (call  D0PRNG81 ) 

2.  a)  Calculate  the  range  at  each  range  gate:  R(I)  =  BRG+I-DRG 
b)  Calculate  log^Q  (R(I)^) 

2.  Calculate  PRT  (pulse  repetition  time) 


Figure  A.ti  Flint'  ,-hurl 
PRELIM. 


t.hc  Eul'ro'u  hi i 
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SUBROUTINE  UNFOLD 


UNFOLD  checks  velocities  at  each  range  gate  after  the  first  to  see  if 
the  velocity  is  aliased,  and  corrects  aliased  velocities  (Figure  A. 7).  The 
process  below  is  repeated  for  the  2nd  to  the  1524^  range  gate. 

1.  ONE:  current  velocity 

OLD:  velocity  at  preceding  range  gate 

2.  If  ONE  =  999  or  OLD  =  999  or  R(I)<0,  go  to  6 

3.  DIFONE  =  | ONE-OLD [ 

4.  If  ONE-’O  VEL=0NE-64* 

If  0NE<0  VEL=0NG+64* 

If  0NE=0  VEL=0 

DIFTW0= ! VEL-OLD j 

5.  If  DIF0NE>DIFTW0  current  velocity  =  VEL 

6.  Stop 


Figure  A. 7  Flew  chart  of  the  Subroutine 
UNFOLD. 


*At  this  point  in  the  program  VEL  can 
have  64  digital  categories  that  span  the 
Nyquist  interval.  Depending  on  the 
actual  PRT,  these  categories  are  later 
scaled  into  appropriate  velocities. 
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